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body volume having at least one target volume and at least one non-target volume. The objective function has a first term related to 
the at least one target volume and a second term related to the at least one non-target volume. The optimization further involves 
determining a minimum of the objective function, whereby beams mapped so as to pass through the at least one non-target 
volume are limited such that the second term is zero only if the weights of beamlets passing through the at least one non-target 
volume are zero. This limit helps to avoid the occurrence of negative beam weights, thereby facilitating computationally efficient 
determination of the minimum of the objective function using matrix inversion. Following the optimization, radiotherapy is delivered 
based on the determined minimum of the objective function. 
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(57) Abstract: The invention relates 
to improved methods and systems for 
computationally efficient optimization 
of radiation dose delivery. The 
optimization involves determining an 
improved form of objective function 
to be used for mapping radiotherapy 
beams to a patient body volume having 
at least one target volume and at least 
one non-target volume. The objective 
function has a first term related to 
the at least one target volume and a 
second term related to the at least one 
non-target volume. The optimization 
further involves determining a 
minimum of the objective function, 
whereby beams mapped so as to pass 
through the at least one non-target 
volume are limited such that the 
second term is zero only if the weights 
of beamlets passing through the at least 
one non-target volume are zero. This 
limit helps to avoid the occurrence 
of negative beam weights, thereby 
facilitating computationally efficient 
determination of the minimum of 
the objective function using matrix 
inversion. Following the optimization, 
radiotherapy is delivered based on the 
determined minimum of the objective 
function. 
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Title: METHOD AND SYSTEM FOR OPTIMIZING DOSE DELIVERY OF 
RADIATION 

Field of the Invention 

[0001] The present invention relates to methods and systems for 
optimizing dose delivery of radiation. In particular, the present invention 
relates to efficient and effective methods of determining a minimum of an 
5 objective function for planning dose delivery of radiation. 

Background of the invention 

[0002] For cancer patients, radiation therapy is recognized as a 
valuable form of treatment. Radiation therapy involves the transmission of 
radiation energy to a tumor site within the patient. 

10 [0003] Radiation therapy planning may be carried out according to a 
forward planning technique or an inverse planning technique. Forward 
planning involves delivering an initial planned radiation dose and then 
delivering subsequent doses by observation or inference of the efficacy of the 
preceding dose in a trial-and-error manner. The optimization of dose delivery 

15 by forward planning is therefore performed according to human observation 
and experience. Inverse planning instead seeks to calculate an optimized 
dose delivery and then work backwards to determine the appropriate radiation 
beam characteristics to deliver that optimized dose. 

[0004] Inverse planning of radiation therapy for tumors may be 
20 performed for Tomotherapy or Intensity Modulated Radiation Therapy (IMRT) 
radiation delivery techniques. Both of these techniques involve transmission 
of radiation beams, usually collimated by a multi-leaf collimator (MLC), toward 
the tumor site from various angular orientations. For Tomotherapy, a helical 
arc is employed to irradiate the tumor slice by slice, while for IMRT multiple 
25 intensity-modulated conical beams are used to irradiate the tumor from a 
number of different directions. 

[0005] In order to ensure that the patient is optimally treated, it is 
necessary to ensure that the radiation dose is deposited primarily within the 
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tumor volume, rather than in the surrounding tissue or organs. It has been 
found to be problematic to quickly and reliably determine an optimization so 
as to maximize the dose delivery to the tumor site while minimizing radiation 
dose delivery to other organs or tissues. 

5 [0006] A fast optimization algorithm is important, not only for designing 
good radiation treatment plans, but also for the successful implementation of 
future interactive adaptive treatment techniques. Conventional optimization 
algorithms using numerical searches, such as the known conjugate gradient 
search with positive beam weight constraints, usually require many iterations 
10 involving long computation times and may result in sub-optimal plans due to 
trapping in local minima of the objective function. 

[0007] It is possible to determine a direct solution of the inverse 
problem using conventional quadratic objective functions, without imposing 
positive beam weight constraints. This solution is computationally faster but 

15 results in unrealistic (negative) beam intensities. Once an ad-hoc condition 
requiring the beam intensities to be positive is introduced (i.e., by forcing 
negative intensity values to be zero), the solution of the quadratic objective 
function by linear algebraic equations yields a radiation therapy dose 
distribution with significant artifacts. These artifacts may significantly 

20 deteriorate an otherwise optimized dose delivery. Accordingly, rather than 
treat a patient with a sub-optimal dose delivery, the rather more 
computationally intensive numerical searching has been preferred for finding 
the minimum of the objective function. 

[0008] A further drawback of current IMRT plan optimization, is that 
25 only about seven to eleven different gantry angles may be employed because 
present techniques find it too computationally intensive to optimize the 
objective function for a greater number of beams. 

[0009] In view of the above shortcomings of existing systems, it is 
desired to provide a method and system for optimized dose delivery, which 
30 addresses or ameliorates one or more of the mentioned shortcomings. 
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Summarv of the invention 

[0010] In one aspect, the invention provides a method of dose delivery 
of radiation. The method comprises the step of determining an objective 
function to be used for mapping radiation beams to a body volume comprising 
5 at least one target volume, and at least one non-target volume. The objective 
function comprises a first term related to the at least one target volume and a 
second term related to the at least one non-target volume. The method 
further comprises determining a minimum of the objective function whereby 
beams mapped so as to pass through the at least one non-target volume are 
10 limited such that the second term is zero only if the weights of beamlets 
passing through the at least one non-target volume are zero. Radiation is 
delivered based on the determined minimum of the objective function. 

[0011] In another aspect, the invention provides a method of 
15 determining an objective function to be used for mapping radiation beams for 
a body volume comprising at least one target volume and at least one non- 
target volume. The objective function comprises a first term related to the at 
least one target volume and a second term related to the at least one non- 
target volume. The method comprises determining a minimum of the 
20 objective function whereby beams mapped so as to pass through the at least 
one non-target volume are limited such that the second term is zero only if 
intensities of beamlets passing through the at least one non-target volume are 
zero. 

[0012] In a further aspect, the invention provides a method of providing 
25 radiation, comprising determining an objective function for optimizing radiation 
dose delivery to a target volume, the objective function having a symmetry 
term for enabling symmetrical dose delivery about an axis of the target 
volume, and providing radiation based on the objective function. 

[0013] In yet a further aspect, the present invention provides a system 
30 for optimizing dose delivery of radiation. The system comprises computer 
processing means for determining an objective function to be used for 
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mapping radiation beams to a body volume comprising at least one target 
volume, and at least one non-target volume, the objective function comprising 
a first term related to the at least one target volume and a second term related 
to the at least one non-target volume. The computer processing means is 
5 arranged to determine a minimum of the objective function whereby beams 
mapped so as to pass through the at least one non-target volume are limited 
such that the second term is zero only if the weights of beamlets passing 
through the at least one non-target volume are zero. The system further 
comprises data communication means operably associated with the computer 
10 processing means for providing data to a radiation delivery apparatus for 
delivering radiation to the body volume based on the determined minimum of 
the objective function. 

[0014] In a still further aspect, the invention provides computer 
readable storage having stored thereon computer program instructions 

1 5 executable on a computer system for causing the computer system to perform 
a dose optimization method. The dose optimization method comprises 
determining an objective function to be used for mapping radiation beams for 
a body volume comprising at least one target volume and at least one non- 
target volume, the objective function comprising a first term related to the at 

20 least one target volume and a second term related to the at least one non- 
target volume; and determining a minimum of the objective function whereby 
beams mapped so as to pass through at least one non-target volume are 
limited such that the second term is zero only if intensities of beamlets 
passing through the at least one non-target volume are zero. 

25 [0015] In the above aspects, it is preferred that the radiation delivery be 
by IMRT or Tomotherapy. 

[0016] Advantageously, embodiments of the invention enable the 
objective function to be minimized according to the solution of a set of linear 
algebraic equations. While there are a number of ways to solve a system of 
30 linear equations, the preferred method is based on determining the inverse of 
a beamlet intersection matrix. Because of the method of determining the 
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minimum of the objective function described herein, determining the inverse of 
the beamlet intersection matrix greatly reduces the possibility of generating 
anomalous negative beam weights for the beamlets. Accordingly, the 
problems associated with negative beam weights and the constraints imposed 
5 on optimization methods to avoid them may be obviated. 

[0017] Further advantageously, because the technique employed by 
the invention allows the optimization to be framed as a solution of algebraic 
linear equations, the lengthy processing time required to search for the global 
minimum of the objective function is substituted with a significantly improved 

10 processing time. This increase in processing efficiency is measurable in 
orders of magnitude. For example, the present technique can accomplish in 
seconds or minutes what would take several hours with some prior art 
techniques. Accordingly, with methods and systems according to the 
invention, medical staff can greatly reduce the time required for radiation 

15 therapy planning while providing a highly optimal dose delivery plan. 

[0018] Further advantageously, embodiments of the invention enable a 
larger number of radiation delivery angles to be employed, compared with 
previous IMRT techniques. This is due to the high computational efficiency 
with which the optimization is carried out according to the invention, providing 
20 higher quality conformal dose distributions to the tumor site and better quality 
optimizations in avoiding radiation delivery to organs at risk and other organs 
or tissues not forming part of the target site. 

[0019] These and further features of aspects and embodiments of the 
invention will be described in detail in the following detailed description. 

25 Brief description of the drawings 

[0020] Figure 1 is a flow chart of an overall process used in dose 
delivery optimization according to an embodiment of the invention; 

[0021] Figure 2 is a flow chart of an organ and contour digitization sub- 
process of the overall process shown in Figure 1; 
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[0022] Figure 3 is a flow chart of a beams and beamlets set-up sub- 
process of the process shown in Figure 1; 

[0023] Figure 4 is a flow chart of a dose calculations sub-process of the 
process shown in Figure 1; 
5 [0024] Figure 5 is a flow chart of a pre-optimization calculations sub- 
process of the process shown in Figure 1; 

[0025] Figure 6 is a flow chart of a beamlets optimization sub-process 
of the process shown in Figure 1; 

[0026] Figure 7 is a flow chart of a results and statistics output sub- 
1 0 process of the process shown in Figure 1 ; 

[0027] Figure 8 is a block diagram of a system for dose delivery 
optimization according to an embodiment of the invention; and 
[0028] Figures 9A to 9D show example optimized dose distribution 
maps and corresponding dose-volume histograms. 

15 Detailed description of the invention 

[0029] The present invention generally relates to methods and systems 
for optimizing dose delivery of radiation therapy to tumor sites within a patient. 
Typically, the radiation will be directed toward a single tumor site, although it 
is not uncommon for multiple tumor sites to be treated simultaneously. While 

20 embodiments of the present invention are capable of taking into account 
multiple tumor sites, which are encompassed by the planning target volume 
(PTV), for simplicity of description, embodiments will primarily be described as 
they apply to a single PTV. Similarly, only a single organ at risk (OAR) and a 
single volume of other intervening tissues or organs, which is called herein all- 

25 the-rest (ATR), are described. 

[0030] The number of PTV, OAR and ATR volumes, as well as the size 
and relative orientation thereof, will vary from patient to patient and according 
to the desired treatment plan determined by the radiation oncologist. For 
example, it is possible that the optimization may not have to take account of 
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an organ at risk, or the PTV may be entirely within the organ at risk, with very 
little ATR volume to take account of. 

[0031] It will usually be necessary or at least advisable for the 
supervising radiation oncologist or other suitably qualified medical personnel 
5 to determine one or more importance parameters in the objective function by 
which relative importance may be attributed to certain tissue or organ volumes 
within the patient relative to the other volumes. 

[0032] During radiation dose delivery, radiation scattering commonly 
occurs due to the passage of the radiation through the body volume. In this 
10 application, scattering effects are not addressed in the optimization method 
and system. It will be understood that embodiments of the invention can be 
modified to account for scattering effects. 

[0033] The most fundamental requirements of a radiation treatment 
optimization are: (/) the dose is homogeneously deposited in the PTV; (//) the 

15 dose deposited in any OAR does not exceed a threshold value and ideally 
should be zero; (///) the dose deposited in ATR organs and tissue not included 
in the PTV and OARs, should be as small as possible and ideally zero to 
minimize the risk of secondary carcinogenesis; (/V) the dose gradient crossing 
the PTV boundaries should be as high as possible. 

20 [0034] Optimizations are pursued by the minimization of a positive- 
definite objective function, also sometimes termed a "cost function". A 
successful optimization will yield a global minimum to this objective function in 
a short computation time with physically achievable beamlet intensities (i.e. 
greater than or equal to zero). 

25 [0035] The optimized objective function should minimize the dose 
deposited in the ATR and the OARs. Theoretically, the dose in these should 
be zero, although it can never actually be zero in the ATR. Consider a 
simplified example of two beamlets only, one with weight a and the other with 
weight 6. The terms in the traditional optimization function for the ATR and 

30 the OAR are each of the form: 
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P(a+bf, 

where p is the importance parameter of the term. The optimization searches 
then for the minimum: 

min{(a + &) 2 }, 

5 [0036] The minimum is zero and it occurs for a = -b. In other words, 
one of the weights will be negative. This is the result from the solution of a 
linear system of equations. 

[0037] The minimum corresponds to: 

at 

|(fl+6) 2 = 0 

1 0 which results in a + b = 0. 

[0038] The current approach is then to solve instead: 

min{(o+6) 2 } with the constraints a>0, b>0. 

[0039] This can only be solved through a numerical search. In order to 
address this problem, embodiments of the invention use instead a term of the 
15 form: 

p(a 2 +b 2 ) 

[0040] This term cannot be zero by one beamlet having a negative 
intensity to cancel the other. For this term to be zero, each beamlet intensity 
must be zero. 

20 [0041] The system of equations is obtained from: 

£{<?+b 2 )=o 

dx 

|(« 2 + 6 2 )=0 

that results in 

a=0 



WO 2005/057463 



CA 02587587 2007-05-15 



PCT/CA2004/002108 



-9- 
b=0 

at the minimum. For a quadratic optimization function, there is only one 
minimum, which is the absolute minimum. 

5 [0042] A typical objective function O satisfying the optimization 
conditions stated above is of the form: 

O = PptvOptv + PoarOoar + PatrOatr 

where the pk are importance coefficients (also called importance parameters) 
and the objectivity terms are: 

10 

xePTV I 

all-beamlets 

0 0A R= Z ( Z *>A(X))\ 
xeOAR I 

15 

all-beamlets 

and 0 Am = Z ( Z M(x) ) 2 , 

where w, is the weight of beamlet /, d,- is the dose deposited at destination 
20 point x by beamlet ;' and (f™ is the dose prescribed to the PTV. 

[0043] The main reason for the appearance of negative weights upon 
optimization of the objective function O is the fact that it is usual to require the 
satisfaction of two conflicting demands: on one hand it is required that 0 A tr = 
0 and on the other hand it is necessary for radiation to pass through the ATR 
25 (and possibly OARs) to reach the PTV. A better requirement is that Oatr 
should be zero only if the weights of all the beamlets passing through the ATR 
are zero, as described in the simplified example above. This requirement is 
satisfied if instead of Oatr we use a new ATR term of the form: 

0 Am = Z a "'T ! ^fdf(x). 

xeATR I 
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Similarly for the OAR: 

_ all-beamlets 

5 O aAR = D E wfdf(x). 

xeOAli f 

[0044] The PTV term in the objective function O cannot be written in 
this way. Accordingly, the medical personnel performing the optimization 
needs to set an importance parameter large enough on the ATR and OAR 
10 terms to balance the PTV term. Even small values of importance parameters 
for the OAR and ATR have been found to be sufficient. 

[0045] In a preferred embodiment, a term is added to the objective 
function O that replaces the unrealistic zero-limit for the beamlet weights with 
an equal-weight limit (which will be referred to herein as circular symmetry), 
15 which is usually the initial set of weights before optimization. This term can 
assume different forms, as part of a general family of symmetry terms. This 
term can be in one of the forms: 



or 

20 or Osym 
or 

or Osym 

or Osym 

or other forms satisfying the condition in the next section. In the following we 
25 will use, for illustration, the first form of Osym. 

[0046] With the weights normalized to: 



Osym = (w r 1) z 

all-begmlels 

Osym = Z (w/ 2 -w/) 

all-beamlets 
all-beams inside beam i 

= E I ( I W;)-1] 2 



Osym 



all-beamlets 

= I Z wfdfo) 

xzall-contours I 

all-beamlets 

= Z Z w, 2 di 2 (x) 
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all-bamilets 

2u w, - total number of beamlets, 

/ 

Osym is positive and its minimum is zero when w,- = 1 for all /. 
[0047] The objective function O then becomes: 

5 0=p PTV 0 P7Y +p 0AR OoM+p Am Oatr+ PsruOsyu 

[0048] The underlying approach behind current optimization techniques 
is to start from zero weights and analyze the results as the weights of each 
beamlet are increased. As a result, searches for a minimum do not 
necessarily result in symmetric dose depositions, even when the system 

10 treated may have a symmetry (eg. symmetric under a reflection). The 
symmetry term introduced here, in essence, starts the analysis of the weights 
from the opposite end: with all beams having the same weight. Given that the 
radiation source travels around the isocentre (i.e. the designated centroid of 
the form or volume) describing a circle on each slice, this requirement starts 

1 5 the analysis from a circularly symmetric perspective. 

[0049] The circular symmetry term Osym has been found to introduce a 
high degree of stability in the results, even when coupled with a small 
importance parameter Psym- Moreover, it tends to smooth the dose 
distribution within the body volume, avoiding local hot or cold spots. 

20 [0050] This introduced symmetry term provides a significant bias 
against generation of negative beamlet intensities during minimization of the 
objective function using matrix inversion. This can be observed if, for 
example, all importance parameters apart from Psym are zero. In such a case, 
the optimization of the objective function would yield a plan where all 

25 beamlets have the same unit weight. Thus, a non-zero value for psym biases 
the beamlet weights towards a unit weight distribution. This bias is small for 
small values of psym and is stronger for larger values. If one were to 
iteratively test and observe the beamlet weight distribution, starting with a 
large value of psym and decreasing it in steps, the distribution would resolve 
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from one in which all weights are substantially the same to a distribution in 
which the beamlet weights are substantially optimized, while keeping all 
beamlet weights positive. 



5 function is that, for a contour having a point or axial symmetry around the 
isocentre, the beamlet weight distribution (and hence dose deposit) as a 
function of gantry angles, will closely follow that symmetry. This ability to 
follow symmetries is derived in part from the large number of gantry angles 
which can be accommodated in the optimization method described herein and 
10 translates into an ability to provide high quality conformal dose deposit 
mapping for target volumes in general. 

[0052] With the symmetry term included, a simplified form of the new 
objective function can be expressed as: 



[0051] 



An advantageous effect of the symmetry term in the objective 




15 where 




where d° rsa " k is the dose prescribed to organ k, and 




20 [0053] The optimal set of weights is obtained by minimizing the 
objective function. The minimum occurs when 



— =0 for all wi 

3v j ' 
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organs organs 

with without 

reqired required „ 

so ^ ^d&r no-^ eor^ , n <xr* 



aldose allbeamlelsl | __ 

where ^- = 2 &J £ 4(*X(*) h 2 ^"' I rf/x), 

i \xeorgan k ) xeorgan k 

5 ^^ =2w (.s, <, ' w )" 2 "''''tR,s.'' ,(ihw ) >< 

an, ^L.,(-,-l)..[(^^)-l] 

where is a unit matrix (i.e., a square array with all elements zero except for 
the diagonal elements that are all equal to one). 
1 0 [0054] Calling now 

Z 4(*H(*)and 

xeorgan t 

por S an t ^ d or S an k £ ^ ^ fo| . each beamIet across th g 

xeor&m k 

whole body volume: 

a// all 

mr as 

required required 



15 and pi ose pr gank +p s 
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[0055] With the noted modifications to the objective function, the 
optimization problem for all the beamlet intensities is reduced to the solution 
of a linear system of equations of the form: 




5 where Wj is the (unknown) weight or intensity of beamlet j, ft is a vector 
(referenced herein as the beamlet dose deposit vector or array) of coefficients 
that depends on the dose deposited by beamlet / within the PTV, and ay is a 
matrix (referenced herein as the beamlet intersection matrix) that describes 
the product of the doses deposited by the intersecting pairs of beamlets /' and j 
10 on all organs. 

[0056] The set of optimal beamlet weights is obtained, for example, 
from (1) by inversion: 

allbeamlets 

j 

[0057] Thus, the solution to the (large) system of linear equations (1) 
15 can be obtained quickly and accurately by inverting the matrix a v using 
standard matrix inversion routines and summing the product of inverted matrix 
ay 1 with vector $for each beamlet /. 

[0058] In one embodiment, the importance parameters for each region 
can be generalized to be region-dependent, i.e. to have different values in a 
20 region within a contour, in which case they can be of the form 

Pregion k ~ Pregionflregionk^ 

where p is an overall constant and q(x) is a function of position. 

[0059] In this case, the definitions of the arrays a and p may be 
generalized by 



25 



xeregionk 
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and fi**=cr*» £ q^Wdjix) 

xeregiotik 

Using these definitions, the matrix formulation of the optimization process now 
becomes: 

a tJ = Pprv^ + (poAR^ + PATR<xf R + Psym ) 
5 fij = J) PTV fij™ + p^ . 

[0060] The previous linear system of equations is now 

and the optimized solution is obtained by the inversion: 

j 

10 [0061] If the functional dependence of the importance coefficients 
remains unaltered, a search of the best set of importance parameters is 
reduced to a search of the best set p regionk in which case the arrays d regionh 

and P regionk do not need to be recalculated. Matrices a i} and are obtained 
from the matrices and by simply setting q regionk = 1 for all regions, 
15 returning to the case in which the importance coefficient has a unique value 
within a contour. 

[0062] Turning now to the drawings, Figure 1 is a block diagram 
illustrating an optimization process 100 according to an embodiment of the 
invention. The optimization process 100 involves obtaining scanned input 

20 data 110 from a scanning apparatus which outputs a series of scans, for 
example such as computed tomography (CT) scans. This scanned input data 
includes a series of "slices" through the body. Each of these slices shows a 
part of the tumor volume in cross-section, together with the remaining body 
volumes, including any organs at risk. When these slices are aggregated as a 

25 series of parallel slices, a three-dimensional picture of the target tumor 
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volume and other body volumes can be obtained. Accordingly, the input data 
110 includes data concerning a number of such parallel slices, sufficient to 
describe the body volume, including the PTV, to which radiation will be 
directed. The input data 110 may be in the Dicom RT standard file format 
5 (including standardized radiation therapy-specific data objects), which can be 
generated by most CT scanning systems or other treatment planning 
systems. Further details of the Dicom RT standard can be obtained from the 
National Electrical Manufacturers Association (NEMA). 

[0063] The input data 110 is received by an optimization module 120, 
1 0 which processes this input data, as described further in relation to Figures 2 to 
7. Once the optimization module 120 has processed the input data 110, an 
output file 130 of the optimized beamlet intensities is generated and output to 
a radiation dose delivery apparatus, such as a medical linear accelerator, 
through a suitable Dicom RT protocol. The radiation therapy can then be 
1 5 delivered according to the optimized dose delivery. 

[0064] Output file 130 is formatted so as to provide sequencing data for 
mapping the optimized beamlet intensities to the leafs of a multi-leaf 
collimator. This is done according to existing leaf sequencing algorithms. 

[0065] The scan data 110 can be stored (eg., in memory 20, shown in 
20 Figure 8) and used to perform several dose optimizations over a period of 
time, as the optimization process facilitates adaptive adjustment of dose 
delivery planning based on different user input requirements. 

[0066] Optimization module 120 is comprised of a series of computer 
program instructions aggregated to form a software program executable by a 
25 computer system (such as computer system 12 in Figure 8, described later). 
Optimization module 120 is configured to receive the scan data input 110 in a 
known file format (e.g., Dicom RT) and to provide the optimized beamlet 
intensities in output file 130 in a corresponding known file format. 

[0067] Optimization module 120 executes a number of sequential 

30 steps, grouped as several sets of steps, which are referred to as 
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subprocesses, as part of the overall optimization process 100. These 
subprocesses include organ and contour digitization 200, beam and beamlet 
setup 300, dose calculation 400, pre-optimization calculation 500, beamlet 
optimization 600 and results and statistics generation 700. These 
5 subprocesses are described in further detail below, with reference to Figures 
2 to 7. 

[0068] Referring now to Figure 2, organ and contour digitization 
subprocess 200 is described in further detail. Organ and contour digitization 
subprocess 200 handles the input of contour data from the CT scans and 

10 digitizes the contours so as to map them onto an underlying grid of cells to 
which all parts of the body volume shown in the CT scan are mapped. 
Subprocess 200 establishes a main array for storing data for all of the cells in 
the body for each CT slice. Subsets of the elements of the main array are 
also stored, corresponding to all of the contours within the body, such as the 

15 PTVandOAR. 

[0069] At step 205, the input data 110 is received as input to the 
optimization module 120. The input data 110 received in this respect includes 
information including, for example, the treatment type (eg. Tomotherapy, 
IMRT) to be performed, the required dose to be deposited within the PTV, the 

20 dose-volume constraints and the CT scans, including organ contours, 
determined by the radiation oncologist. The dose-volume constraints indicate 
the maximum or minimum radiation dose to be delivered to a particular 
volume. For example, for an organ at risk such as the spinal cord, a 
constraint may be provided such that no more than 45Gy of radiation dose 

25 should be received in any part of the OAR volume. 

[0070] Once the input data 110 is received at step 205, supervising 
medical personnel may choose, at step 210, the grid cell size and resolution 
to be used for optimizing the radiation dose delivery. This information is then 
used to generate a discretized grid having cells of the chosen size. For each 

30 CT scan, all organs, contours and beams are mapped onto a single main 
array representing grid cells of the chosen spatial size within the body volume. 
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Arrays representing cells with each organ contour, beam and beamlet are 
subsets of this main array. Each element of each array includes the grid 
coordinates of the corresponding cell in the grid. A typical cell size employed 
by embodiments of the invention may be 1mm square in the plane of each 
5 slice. 

[0071] The properties of an organ are assumed to be uniform within 
each cell. Within the descretized grid, beam propagation is calculated with an 
accuracy given by the resolution, which is usually about 1.25 times the width 
(which equals the height) of each cell. 

1 0 [0072] Importantly, the resolution is set greater than the cell size so that 
each beamlet always traverses at least the center of one cell in the grid at 
each depth level of its propagation. This condition produces beamlets that 
are more regular in shape and avoids the beamlets being discontinuous with 
adjacent beamlets. Within each step and within a width equal to the 

15 resolution, the beam properties are assumed to be uniform. 

[0073] Because of the higher computational efficiency enabled by 
embodiments of the invention, more data can be handled by the optimization 
process 100 and a relatively high resolution and small cell size can be 
achieved for the cells of the PTV, OAR and ATR volumes, leading to a more 
20 optimized treatment plan for the patient. 

[0074] At step 215, a contour of the PTV is retrieved from the Dicom 
RT input data 110. At step 220, the contour is checked for clockwise or 
counter-clockwise orientation. If the points of the contour are in a counter- 
clockwise order, the order of those points is reversed so as to be clockwise at 
25 step 225. If the points of the contour are in a clockwise order, the contour is 
digitized so as to interpolate a continuous contour outline from contour points 
provided by the radiation oncologist and the cells within the contour are 
tagged and saved as such, at step 230. 

[0075] The method for determining the orientation of the contour is as 
30 follows. Each slice of each organ is represented by a two-dimensional 
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contour in the plane of that slice. Here we refer to the two-dimensional 
contours in a specific slice. For each contour, the input data specifies a set of 
points or vertices (e.g. x-y coordinates) that outline that contour. These 
vertices are generated by the radiation oncologist on the basis of CT scan 
images. In order to be able to find which points are inside or outside that 
contour, it is necessary to first find out if the set of vertices follows a clockwise 
or anticlockwise direction. 

[0076] Assuming a set of orthogonal axes defined in the plane of the 
slice: 

1 ) Find vertex A: the topmost vertex in the contour. 

2) Find vertex C: the bottom vertex in the contour. 

3) Find vertex B: the rightmost vertex in the contour that is neither 
AorC. 

4) Find vertex D: the leftmost vertex in the contour that is neither A 
orC. 

[0077] Given that the contour encloses a finite area, at least three of 
the above vertices must be distinct. If the contour is in the clockwise direction, 
then any three distinct vertices of the above must be in the order A-B-C-D (or 
any cycle of this order, such as D-A-B-C). If it is not, then the contour is 
determined to be in an anticlockwise orientation and the order of the elements 
of the contour array is inverted to assume clockwise order. 
[0078] After step 225 or 220, digitizing is performed at step 230, 
starting from the topmost vertex and proceeding clockwise to join each 
consecutive pair of vertices by lines. The area inside the contour is effectively 
divided into horizontal lines, each starting at the left boundary and ending at 
the right boundary. As these lines are drawn between vertices, the lines are 
digitized into cells or "pixels" which are inserted in the main array into which 
all contours are digitized. 

[0079] In broad terms: 
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1) If the line is drawn in a direction going downwards, then 
each pixel on that line is a right boundary of the horizontal line at 
the height of the pixel. 

2) If the line is drawn in a direction going upwards, then 
5 each pixel on that line is a left boundary of the horizontal line at 

the height of the pixel. 

3) The position of all left and right boundaries for each line 
at each height within a contour is stored in memory. 

4) After the boundary has been completely digitized, all of 
10 the cells are labeled in each horizontal line between the saved 

left and right boundaries as belonging to the surface enclosed 

by the relevant contour. 
[0080] Steps 215 to 230 are repeated for each contour and for each 
organ volume (eg. OAR, ATR, PTV). Once the last contour has been 
15 digitized, organ and contour digitization subprocess 200 feeds into beam and 
beamlet setup subprocess 300 at step 305, as indicated by reference 
indicator 1 in Figures 2 and 3. 

[0081] Referring now to Figure 3, beam and beamlet setup subprocess 
300 is described. Beam and beamlet setup subprocess 300 determines the 
20 boundaries of each beam as projected from each different angle toward the 
PTV. Subprocess 300 divides the beams into beamlets, calculates the 
beamlet boundaries and determines which cells are inside each beamlet for 
each beam. 

[0082] Subprocess 300 begins with data input from the user (i.e. 

25 medical personnel planning the optimized dose delivery) as to the desired 
radiation beam characteristics, at step 305. For example, for IMRT, multiple 
beams (for example, up to 51) may be specified at various angles relative to 
the PTV isocentre. At this step, the user also inputs beam setup information, 
such as the distance between the radiation source and the PTV isocentre. 

30 For Helical Tomotherapy, step 305 commonly involves choosing 51 beams at 
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equally spaced angles, with beamlets separated by about 6.25mm at the PTV 
isocentre. 

[0083] For radiation therapy, the radiation beams are typically delivered 
to a patient lying on a bed while a gantry carrying a radiation beam emitter 
5 moves around the patient. The gantry can be positioned at numerous different 
angles, depending on the dose delivery plan developed by the radiation 
oncologist and the limitations of the radiation delivery apparatus. For IMRT, 
the chosen gantry angles can be along a circular arc around the patient in a 
single plane or in multiple intersecting planes, treating all tumor slices 
10 simultaneously, one gantry angle at a time. For Tomotherapy, each slice is 
treated from set gantry angles along a predetermined circular or helical arc. 

[0084] Once the beam characteristics are chosen at step 305, the first 
of the predetermined gantry angles is set at step 310. The beam boundaries 
are then calculated at step 315, so as to only coincide with the outer-most 
15 edges of the PTV contour, based on the beam setup information and PTV 
contour data. 

[0085] At each gantry angle, the position and width of the beam is 
calculated in order to fully cover the PTV as seen by the beam source from 
the radiation beam emitter. 

20 [0086] Once the planned beam boundaries are determined at step 315, 
the beam is divided into beamlets at step 320. The number of beamlets within 
each beam will depend on the tumor-shape, gantry angle, equipment 
limitations, beam boundaries and previously specified resolution. Once the 
number of planned beamlets is determined, the positions of each of the leafs 

25 in the multi-leaf collimator are calculated for the beam from the selected 
gantry angle. 

[0087] At step 325, for each beam and beamlet, the cells within each 
organ and contour through which each beam and beamlet would pass are 
stored. This allows quick calculation of beam statistics for each beam, as well 
30 as quick calculation of the optimization arrays for the beamlet intensities. 
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[0088] At step 330, subprocess 300 checks whether the last of the 
predetermined gantry angles has been selected, and if not, the next gantry 
angle is selected at step 335. Steps 315 to 325 are repeated for each gantry 
angle chosen at step 335. 
5 [0089] After the last gantry angle has been processed at step 330, 
subprocess 300 feeds into dose calculation subprocess 400, beginning at 
step 405, as indicated by reference indicator 2 in Figures 3 and 4. 
[0090] Referring now to Figure 4, dose calculation subprocess 400 is 
described. Dose calculation subprocess 400 simulates the propagation of 

10 each beamlet toward the PTV so as to determine the amount of radiation 
energy which would be deposited by each beamlet in each cell of each organ, 
based on the previously determined beam and beamlet setup. The 
determined dose deposit in each cell is then stored for later reference. 
[0091] At step 405, the first gantry angle is selected from which 

15 beamlets are to be propagated. At step 410, equal intensity beamlets are 
propagated (as a simulation for planning calculation purposes only) according 
to the beam and beamlet setup determined in subprocess 300. The equal 
intensity beamlets propagated at step 410 have a default normalized 
weighting of 1 . 

20 [0092] For radiation therapy planning calculations, the way in which 
each beamlet of each beam propagates through tissue and deposits energy in 
each cell is calculated as follows. For the purpose of calculating beamlet 
propagation, each (actual) beamlet is (computationally) divided into narrower 
sub-beamlets, termed here "elementary propagators". The width of each of 

25 these elementary propagators (at the isocentre) is equal to the resolution 
(approximately 1 .25 times the linear cell dimension). 

[0093] The energy deposit during propagation of each elementary 
propagator is calculated (according to a known formula) in small steps in the 
direction of propagation according to the resolution. The elementary 
30 propagator is divided or resolved along its length into small trapezoids (due to 
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divergence of the beamlets from the emitter) of linear dimension equal to the 
resolution (i.e. slightly larger than the linear dimension of the cells in the 
underlying grid). Although the linear dimension (i.e., the distance between the 
parallel sides) of each trapezoid is uniform, each succeeding trapezoid is 
5 slightly wider than the last, so that the trapezoids are non-uniform in size. 

[0094] When the center of a cell (in the main underlying grid) lies inside 
one of these small trapezoids, the elementary propagator is determined to 
deposit energy in that cell. By making the resolution slightly larger than the 
linear dimension of a cell, as the elementary propagator propagates, every 
10 single trapezoid into which it is divided can be considered to deposit energy 
into at least one cell, thus making the elementary propagator continuous and 
not fragmented. 

[0095] The proportion of the resolution to the cell width may vary, 
depending on requirements, but is preferably between 1 and 2 times the cell 
15 width. 

[0096] At step 415, the dose to be deposited in each cell by the 
propagated beamlets is stored for each cell of each organ or body volume 
affected by the beamlets of the beam at the selected gantry angle. 
[0097] Dose calculation subprocess 400 checks at step 420 whether 
20 the selected gantry angle is the last angle at which simulated beamlet 
intensities are to be propagated and, if it is not, the next gantry angle is 
selected at step 425 and steps 410 and 415 are repeated until the last gantry 
angle has been simulated. 

[0098] The planned beamlet intensities stored as part of subprocess 
25 400 are used to calculate the dose which would be deposited by each 
beamlet in the beam at the selected gantry angle for each cell at each volume 
affected by the beamlets. The dose which would be deposited by each 
beamlet is calculated by propagating several elementary propagators of 
radiation per beamlet (as described above), each propagator being of equal 
30 width (at the isocentre) to the resolution. This calculation is performed 
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according to existing dose deposition formulae. In storage step 415, the 
arrays of the cell-by-cell dose deposit data of all beamlets for all organs are 
stored (eg. in RAM) for later use without recalculation. 
[0099] Dose calculation subprocess 400 feeds into pre-optimization 
5 calculation subprocess 500 in Figure 5. Reference indicator 3 in Figures 4 
and 5 joins the respective flowcharts in this regard. 

[00100] Referring now to Figure 5, pre-optimization calculation 
subprocess 500 is described. Pre-optimization calculation subprocess 500 
calculates and stores the matrices and arrays of dose deposit data in each 
10 cell for each beamlet and for intersecting beamlets, in order to perform the 
optimization calculations. 

[00101] Subprocess 500 begins with step 505, at which dose deposition 

coefficient arrays p organ and beamlet intersection matrices a or9an for all organs 

and contours are calculated and stored. Each matrix a or9an comprises 
1 5 elements that are the products of the doses deposited by intersecting pairs of 

beamlets (i.e. from different gantry angles) across all cells in an organ (i.e. 

PTV, ATR or OAR). Array p or9an is a vector of coefficients of the dose to be 

deposited by each beamlet within the cells of each organ. 

[00102] For each organ, k, the matrix a organ has elements ay sa " 1 labeled 
20 by the indices / and determined according to the expression: 

a omo» t= £ d,(x)dj(x) 

where the summation runs over all points x inside organ k and dtfx) and d/xj 
are the doses deposited at point x by unit-weight beamlets i and j, 
respectively. The terms of this matrix are quickly calculated from the arrays 
25 stored in fast memory (RAM) for the unit weight beamlets in step 41 5. 

[00103] As part of step 505, the arrays and matrices for each beamlet 
propagated within each organ are stored in memory (eg. RAM) for later quick 
retrieval during the optimization calculations. The calculations at step 505 are 
performed only once for the initial set of equal intensity beamlets. 
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[00104] At step 510, default importance coefficients for the objective 
function are set. These default coefficients are set according to previous 
experience with appropriate weighting. These coefficients are used to 

achieve a workable optimization of the objective function O. The medical 
5 physicist or other medical personnel performing the optimization may choose 
the default coefficients and may alter these later as part of beamlet 
optimization subprocess 600. 

[00105] Typical default values for the importance parameters are, for 
example, p PT v = 40; poar = 28; patr = 1; and p S y/v? = 1. If a better conformal 

10 dose deposit within the PTV is needed, then Pptv is increased, for example to 
100 or more. Similarly, the values for Poar and Patr are adjusted to suit 
radiation therapy planning requirements. If the system of contours is such that 
with these parameters or with the modified parameters the conditions on the 
PTV or the OAR are excessively demanding and one or more negative 

15 beamlet weights are detected, then one or both of Patr and p$ym can be 
increased, typically to a value like 2 or 3. 

[00106] After pre-optimization calculation subprocess 500, step 510 
feeds into step 605 in beamlet optimization subprocess 600. Figures 5 and 6 
are linked by reference indicator 4 in this regard. 

20 [00107] Referring now to Figure 6, beamlet optimization subprocess 600 
is described. Beamlet optimization subprocess 600 performs the optimization 
calculations by matrix inversion of the beamlet intersection matrix and 
determines the optimal beamlet weights according to the solution of a system 
of linear equations. Beamlet optimization subprocess 600 further determines 

25 the planned dose distribution among all contours and generates graphical 
outputs of the dose distribution. 

[00108] Subprocess 600 begins with step 605, at which the overall 
optimization matrix a and array B are calculated according to the 
precalculated arrays and matrices for the organs and beamlets, taking into 
30 account the predetermined importance coefficients. 
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[00109] Matrix a and array p are calculated for each beamlet by the 
simple expressions: 

allPTV 

A= I P PTVk Pr k P Tv k +Pv m and 
k 

allPTf' allOAR alMm 

k k k 

where <% is a unit matrix and k is the number of contours of each kind (eg. 
OAR, PTV, ATR). This is a fast calculation for which the arrays saved in fast 
memory in step 505 are used. 

[00110] At step 610, the matrix a is inverted using a known fast matrix 
inversion algorithm, such as the lower-upper-diagonal (LUD) algorithm. Any 
suitably computationally efficient matrix inversion algorithm may be used at 
step 61 0. 

[00111] At step 615, subprocess 600 checks matrix a for near-linear 
dependence. This check is performed to ensure that there are no redundant 
or nearly redundant beamlets. If near-linear dependence is found, this is 
corrected at step 620 using a singular value decomposition (SVD) algorithm to 
appropriately condition matrix a. 

[00112] In the absence of near-linear dependence, the beamlet weights 
are calculated, at step 625, as the product of inverted matrix a with dose 
deposit array p. As part of step 625, the optimized beamlet weights are 
determined according to the solution of a linear system of equations (resulting 
from the product of inverted matrix a" 1 with dose deposit array p). This 
solution can be obtained by solving a system of N linear equations in N 
variables, where N is the number of beamlets. 

[00113] At step 630, beamlet optimization subprocess 600 checks 
whether any of the beamlet weights have been calculated at step 625 to be 
negative. If there are any negative weights (i.e. negative beamlet intensities), 
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the user is notified and, at step 635, is advised to adjust one or more of the 
importance coefficients in the objective function. Following adjustment of the 
importance coefficients, steps 605 to 625 are repeated until no negative 
beamlet weights are output from the calculations of step 625. 

5 [00114] If the calculated beamlet weights are positive or zero, the 
optimized dose to be deposited in all organs and contours is recalculated with 
the optimized beamlet weights at step 640. At step 645, dose deposit 
statistics are calculated for all organs, contours and beamlets, including dose 
volume histogram (DVH) plots, for the optimized dose delivery plan. 

10 [00115] Graphics, such as colour-coded dose distribution maps, are 
generated at step 650 according to the calculated dose statistics, where the 
colour-coded dose distribution is overlaid on the contours to provide an easy 
indication of the dose distribution across all contoured volumes. Example 
dose distribution maps and dose-volume histograms are shown in Figures 9A 

15 to 9D. Each colour-coded dose distribution graph is accompanied by a 
corresponding dose volume histogram to provide the user with a more 
accurate indication of whether the dose-volume constraints will be met by the 
proposed optimization of beamlet weights. 

[001 1 6] At step 655, the user is given the opportunity to indicate whether 
20 the planned dose distribution is satisfactory and, if not, is prompted to change 
the importance parameters at step 660. If the user elects to change the 
importance coefficients at step 660, steps 605 to 655 are re-executed until a 
satisfactory dose distribution is achieved. 

[00117] Alternatively, if the user wishes to change some of the physical 
25 setup characteristics, such as the gantry angles, optimization process 100 
returns (not shown) the user to subprocess 200 or 300, as appropriate. 
[00118] Once the dose distribution is considered to be satisfactory, 
beamlet optimization subprocess 600 feeds into step 705 of result and 
statistics generation subprocess 700. Figures 6 and 7 are joined by reference 
30 indicator 5 in this respect. 
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[00119] Referring now to Figure 7, results and statistics generation 
subprocess 700 is described. Results and statistics generation subprocess 
700 generates an output of the optimized beamlet weights after the planned 
dose distribution is approved or refined. 

5 [00120] At step 705, it is determined whether the desired form of 
radiation therapy is IMRT, and if so, generates leaf sequencing data for a 
multi-leaf collimator (MLC), at step 710. If the desired form of radiation 
therapy is Tomotherapy, a beam profile of collimated beamlets is generated 
(not shown) for each beam at each gantry angle and for each tumor slice. 

10 [00121] At step 715, the user is again given the opportunity to refine the 
dose distribution, for example in order to suit the MLC leaf sequencing (if 
IMRT is used) or to accommodate other physical constraints imposed by the 
radiation therapy delivery system. If the dose distribution requires refining, a 
more accurate dose deposition may be substituted for that previously defined 

1 5 and the dose distribution is recalculated at step 720. 

[00122] If no further refinement of the dose distribution is required, an 
output data file is created at step 725, including optimized beamlet intensities 
130 and leaf positions and sequences generated at step 710. Any beam, 
organ or contour statistics, together with data for displaying colour-coded 

20 dose distributions and dose-volume histograms may also be output (eg. to a 
display) at step 725, if desired. The statistics and data for generating dose 
distribution graphs and histograms are stored in a memory of the computer 
system running optimization module 100 for user review and display and the 
optimized beamlet weights and MLC leaf positions and sequences are output 

25 to the radiation delivery system to begin radiation therapy treatment of the 
patient. 

[00123] While the methods and subprocesses for optimization have 
been described above in relation to one preferred embodiment, the invention 
may be embodied also in an optimization system 10, running optimization 
30 module 120, as is shown in Figure 8. 
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[00124] In Figure 8, optimization system 10 includes a computer system 
12 having memory 20 and optimization module 120 running as executable 
computer program instructions thereon. Memory 20 is fast memory, such as 
fast-access RAM, for storing arrays and matrices and calculation terms used 
5 during the optimization process 100. Memory 20 is also used to store 
statistics and/or calculations for generating dose-volume histograms and 
colour-coded dose distribution graphics. Memory 20 may include, or 
communicate with, secondary (slower) memory to facilitate appropriate data 
storage during process 100. Optimization module 120 uses memory 20 as 
1 0 required for its storage requirements. 

[00125] Computer system 12 further includes normal computer 
peripherals, including graphics displays, keyboard, secondary memory and 
serial interface, as would normally be used for a computer system which 
receives input data 110 and generates corresponding output data 130. 

15 [00126] In a further aspect, the invention may be embodied in computer 
program instructions (i.e. software for executing the described methods) 
stored on computer-readable storage media, such as a hard disk, CD-ROM or 
RAM. 

[00127] While preferred embodiments of the invention have been 
20 described in relation to dose delivery of radiation for radiation therapy 
treatment, it is to be understood that the optimization process 100 and 
optimization system 10 may be equally useful for planning optimized radiation 
delivery to body volumes other than those of human patients under treatment 
for cancerous tumors. For example, the described systems and methods may 
25 be employed for animals other than humans and may be employed for 
irradiating non-living tissue or material or organic matter where selective dose 
delivery of radiation is desired. 

[00128] Figures 9A to 9D show example optimized dose distribution 
maps 902 and dose-volume histograms 904 for IMRT. Figures 9A and 9B 
30 show a dose distribution map 902 of a CT scan slice of a neck tumor. Dose 
distribution map 902 shows an outer body contour 913 enclosing all volumes 
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of interest. OAR contour 915 encloses an OAR volume 914, which in this 
example is the spine. PTV contour 917 encloses a PTV volume 916, which is 
the neck tumor. The remaining volume within outer body contour 913, which 
is not within OAR and PTV volumes 914 and 916, is called the ATR volume 
5 912. 

[00129] As the example in Figure 9A illustrates, it can be difficult to 
direct beams so as to deliver radiation to the PTV volume 916 without also 
directing some radiation towards the OAR volume 914. In this example, it is 
particularly important to minimize delivery of radiation to the spine as it is 

10 sensitive to radiation delivery and an excessive dose may result in damage to 
the spinal cord or nerve endings therein. Typically, dose-volume constraints 
for an OAR such as the spinal cord are such that none of the OAR volume 
should receive a dose in excess of about 45 Gy. If only a low number of 
gantry angles are employed in the optimization planning, this dose-volume 

15 constraint may not be able to be met, whereas if a larger number of gantry 
angles are employed, such a dose-volume constraint can be met. The greater 
computational efficiency achievable by the invention enables a larger number 
of gantry angles to be employed, which results in a better conformal mapping 
of the dose delivery plan. 

20 [001 30] Figures 9A and 9B show the same ATR, OAR and PTV volumes 
and contours but, whereas Figure 9A shows the output plan of an optimization 
using twenty gantry angles, Figure 9B shows a plan using sixty gantry angles. 
In Figure 9B, the OAR curve 924 of dose-volume histogram 904 shows that 
none of the OAR volume 914 receives more than 20% of the dose, as 

25 compared to the 40% indicated in Figure 9A. The increase in the number of 
gantry angles allows for greater flexibility in optimizing the dose delivery plan 
so as to avoid irradiating the OAR while maximizing the radiation dose to the 
PTV. This increase in the number of gantry angles is enabled by the 
increased computational efficiency of the present optimization method. 
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[00131] In dose-volume histogram 904, ATR curve 922 indicates the 
dose-volume distribution of ATR volume 912, while PTV curve 926 indicates 
the dose-volume distribution to PTV volume 916. 

[00132] While the optimization is faster for fewer gantry angles (and 
5 thus fewer beamlets), for example in the order of 5 seconds, it is less 
optimized than the plan using sixty gantry angles in the sense that it delivers 
a higher average dose to the OAR volume 914. Conversely, while the 
optimization using sixty gantry angles is more accurate in avoiding the OAR 
volume 914, the time required for the optimization is greater, for example in 
1 0 the order of 2 minutes and 1 2 seconds. 

[00133] The reference indicators mentioned above in relation to Figures 
9A and 9B are also applicable to Figures 9C and 9D. 

[00134] Figures 9C and 9D show further examples of dose distribution 
maps 902 and corresponding dose-volume histograms 904. Outer body 

15 contour 913 is not shown in Figures 9C and 9D because in this example dose 
distribution map 902 is a close-up within a larger body volume. Figures 9C 
and 9D relate to a planned dose delivery to a prostate tumor, indicated by 
PTV volume 916. In this example, the OAR volume 914 is the bladder. As 
can be seen from Figures 9C and 9D, the OAR body contour 915 overlaps the 

20 PTV contour 917 and therefore it is not possible to minimize dose delivery to 
the OAR to the same degree as would be possible where the contours and 
volumes did not overlap. 

[00135] The optimizations shown in Figures 9C and 9D use the same 
number of IMRT gantry angles and beamlets, with the same beamlet width of 

25 2mm. The primary difference between Figures 9C and 9D is that the 
optimization shown in Figure 9D was designed to provide a highly conformal 
dose distribution to the PTV, with less importance being given to irradiating 
the OAR. A comparison of OAR curves 924 in Figures 9C and 9D indicates 
that the OAR volume 914 in Figure 9D received a higher average dose than 

30 that of Figure 9C, primarily because of the lower importance attributed to 
minimizing irradiation of the OAR in the optimization shown in Figure 9D. 
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Thus, manipulation of importance coefficients can have a significant effect on 
the resultant dose distribution within the PTV and OAR. 

[00136] In this description, certain terms have been used 
interchangeably. For example, beamlet weights and beamlet intensities have 
5 been used interchangeably and are intended to have the same meaning. 
Similarly, importance parameters and importance coefficients have been used 
interchangeably and are intended to have the same meaning. Also, some 
terms used in this description may be called by other names in related 
technical papers, although the meaning is the same. For example, the term 
10 objective function used herein may be called a cost function by others. 
Similarly, the term fluence is used elsewhere for describing beamlet weights 
or intensities. It is intended that this description should cover all terms having 
the same meaning, as would be understood by the skilled person, as the 
terms used herein. Further reference can be made to the Glossary below. 
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Glossarv 

[00137] Beam: a ray of radiation transmitted in two or three dimensions. 

[00138] Beamlet: A fraction of a beam of radiation. A beam can be 
divided into many small beamlets. For example, a beam with a 10x10 cm 2 
5 cross-section can be divided into 100 beamlets with 1x1 cm 2 cross-section. 

[00139] Beamlet Dose Deposit Array: an array of the planned dose 
deposits by a particular unit weight beamlet in the cells within a particular 
contour, multiplied by the prescribed dose for that contour. 

[00140] Beamlet Intersection Matrix: a matrix of the products of the 
10 separate planned dose deposits in a cell by intersecting beamlets added over 
all cells within a particular contour. The diagonal elements of this matrix 
contain the sum of the squares of the dose deposits of each beamlet over all 
the cells of the contour. 

[00141] Beamlet weight: Intensity level of a beamlet. 

15 [00142] Contour: a two-dimensional outline of a structure in a single CT 
scan, outlining, for example, a target volume or an organ or any other internal 
mass. 

[00143] Importance Parameter (or coefficient): a relative importance 
factor used in the optimization process. An importance parameter can be 

20 defined in the range from 0 (ie. where the corresponding organ is considered 
to be irrelevant to the optimization) to any positive value. The higher the value 
of the importance parameter relative to other importance parameters, the 
greater the importance attributed to the corresponding term of the objective 
function and the greater influence that term will have on the optimization of 

25 the objective function. 

[00144] Intensity: a measure of relative strength of radiation. 
[00145] Non-target volume: the volume outside the target volume but 
within the body volume. The non-target volume should receive as small a 
radiation dose as possible. 
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[00146] Objective function: a mathematical function used for optimizing 
the way radiation is delivered to a body volume. 

[00147] Organ-at-risk: an organ inside a patient that may be adversely 
affected or damaged by radiation during radiation therapy. Generally, 
5 important organs near a tumor are organs at risk since they may receive high 
radiation doses. Typical organs at risk are: the spinal cord, eyes, liver, lung, 
kidney, heart, rectum, bladder, brain stem, optic nerves, optical chiasma, 
parotid and small intestine. 

[00148] Planned radiation dose deposit: the calculated radiation dose to 
10 be delivered at a point in a patient for a particular plan. 

[00149] Radiation: Energy radiated or transmitted as rays, waves, in the 
form of particles. For this application, we refer to high energy x-rays 
particularly. 

[00150] Target volume: a contoured volume, for example, including a 
15 tumor and surrounding tissue, to receive a certain prescription of radiation 
dose to treat a tumor or other mass within that volume. 
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Claims: 

1 . A method of dose delivery of radiation comprising the steps of: 

determining an objective function to be used for mapping radiation 

beams to a body volume comprising at least one target volume, and at least 
5 one non-target volume, the objective function comprising a first term related to 

the at least one target volume and a second term related to the at least one 

non-target volume; 

determining a minimum of the objective function whereby beams 

mapped so as to pass through the at least one non-target volume are limited 
10 such that the second term is zero only if the weights of beamlets passing 

through the at least one non-target volume are zero; and 

delivering radiation based on the determined minimum of the objective 

function. 

15 2. The method of claim 1, wherein the second term comprises, for all of a 
plurality of non-target volume portions, a non-target volume sum of beamlet 
sums related to respective non-target volume portions, each beamlet sum 
being a sum of the product of the squared weight of the beamlet with the 
squared planned radiation dose deposit at the respective non-target volume 

20 portion. 

3. The method of claim 1, wherein the objective function further 
comprises a third term related to an organ-at-risk (OAR) volume and wherein 
the third term comprises, for all of a plurality of OAR volume portions, an OAR 
25 sum of beamlet sums related to respective OAR volume portions, each 
beamlet sum being a sum of the product of the squared weight of the beamlet 
with the squared planned radiation dose deposit at the respective OAR 
volume portion. 
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4. The method of claim 1, wherein the objective function further 
comprises a symmetry term for enabling symmetrical dose delivery about an 
axis of the at least one target volume. 

5 5. The method of claim 4, wherein the symmetry term is of the form: 



6. The method of claim 4, wherein the symmetry term is positive and its 
minimum is zero when w/ = 1 for all /', where w, is the weight of beamlet /* of a 
plurality of radiation beams. 

15 

7. The method of claim 1, wherein the step of determining a minimum 
includes solving a linear system of equations to determine the weights of the 
beamlets. 

20 8. The method of claim 7, wherein the solution of the linear system of 
equations is generated using matrix inversion of a beamlet intersection matrix 
for each beamlet. 

9. The method of claim 8, wherein the solution of the linear system of 
25 equations is generated by the product of the inverted beamlet intersection 

matrix with a beamlet dose deposit array. 

10. The method of claim 8, wherein the beamlet intersection matrix 
comprises a sum of organ volume matrices respectively corresponding to the 

30 at least one target volume and the at least one non-target volume, each organ 
volume matrix being weighted by a respective importance parameter. 




10 



where Osym is the symmetry term, and 

w; is the weight of beamlet / of a plurality of radiation beams. 
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11. The method of claim 1 0, wherein the beamlet intersection matrix further 
comprises a symmetry term having a symmetry importance parameter for 
weighting the symmetry term. 

5 12. The method of claim 1 , further comprising: 

receiving contour data relating to a two-dimensional contour of the at 
least one target volume or the at least one non-target volume; 

determining from the contour data whether the contour is oriented 
clockwise or anti-clockwise; and 
10 if the contour is determined to be anti-clockwise, changing the order of 

the contour data so that the contour is oriented clockwise. 

13. The method of claim 12, wherein determining whether the contour is 
oriented clockwise or anti-clockwise further comprises: 
1 5 a) determining a topmost vertex of the contour; 

b) determining a lowermost vertex of the contour; 

c) determining a rightmost vertex of the contour that is neither the 
topmost or lowermost vertex; 

d) determining a leftmost vertex of the contour that is neither the 
20 topmost or lowermost vertex; and 

e) determining the contour orientation according to the relative 
clockwise order of the topmost, lowermost, rightmost and leftmost vertices 
with respect to each other. 

25 14. The method of claim 12, further comprising: 

extrapolating a continuous contour from the contour data; 
determining all right and left boundaries of the continuous contour; and 
determining a cell of the body volume to be within the continuous 
contour if the cell lies between a facing pair of right and left boundaries. 

30 

15. The method of claim 14, wherein a boundary is determined to be a left 
boundary if the contour data indicates an upwardly extending sequence of 
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contour points and a boundary is determined to be a right boundary if the 
contour data indicates a downwardly extending sequence of contour points. 

16. The method of claim 1, wherein said body volume is virtually divided 
5 into a plurality of cells of a predetermined size and said radiation beams are 

mapped to said body volume such that fractions of the radiation beams are 
dimensioned proportionally to the size of said cells. 

17. The method of claim 16, wherein said fractions are resolved into 
1 0 linearly sequential portions of non-uniform size. 

18. The method of claim 17, wherein a linear dimension of said sequential 
portions is uniform and is 1 to 2 times a width dimension of said cells. 

15 19. The method of claim 18, wherein said linear dimension is about 1.25 
times said width dimension. 

20. The method of claim 1, wherein the dose delivery of radiation 
comprises intensity-modulated radiation therapy. 

20 

21. The method of claim 1, wherein the dose delivery of radiation 
comprises Tomotherapy. 

22. A method of determining an objective function to be used for mapping 
25 radiation beams to a body volume comprising at least one target volume and 

at least one non-target volume, the objective function comprising a first term 
related to the at least one target volume and a second term related to the at 
least one non-target volume, the method comprising: 

determining a minimum of the objective function whereby beams 
30 mapped so as to pass through the at least one non-target volume are limited 
such that the second term is zero only if intensities of beamlets passing 
through the at least one non-target volume are zero. 
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23. The method of claim 22, wherein the second term comprises, for all of 
a plurality of non-target volume portions, a non-target volume sum of beamlet 
sums related to respective non-target volume portions, each beamlet sum 

5 being a sum of the product of the squared weight of the beamlet with the 
squared planned radiation dose deposit at the respective non-target volume 
portion. 

24. The method of claim 22, wherein the objective function further 
10 comprises a third term related to an organ-at-risk (OAR) volume and wherein 

the third term comprises, for all of a plurality of OAR volume portions, an OAR 
sum of beamlet sums related to respective OAR volume portions, each 
beamlet sum being a sum of the product of the squared weight of the beamlet 
with the squared planned radiation dose deposit at the respective OAR 
15 volume portion. 

25. The method of claim 22, wherein the objective function further 
comprises a symmetry term for enabling symmetrical dose delivery about an 
axis of the at least one target volume. 

20 

26. The method of claim 25 wherein the symmetry term is of the form: 




25 where Osym is the symmetry term, and 

w/ is the weight of beamlet / of a plurality of radiation beams. 



27. The method of claim 25, wherein the symmetry term is positive and its 
minimum is zero when w ; - = 1 for all /, where w,- is the weight of beamlet / of a 
30 plurality of radiation beams. 
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28. The method of claim 22, wherein the dose delivery of radiation 
comprises intensity-modulated radiation therapy. 

29. The method of claim 22, wherein the dose delivery of radiation 
5 comprises Tomotherapy. 

30. A method of providing radiation, comprising: 

determining an objective function for optimizing radiation dose delivery 
to a target volume, the objective function having a symmetry term for enabling 
1 0 symmetrical dose delivery about an axis of the target volume; and 
providing radiation based on the objective function. 

31 . The method of claim 30, wherein the symmetry term is of the form: 



where Osym is the symmetry term, and 

w, is the weight of beamlet / of a plurality of radiation beams. 

20 32. The method of claim 30, wherein the symmetry term is positive and its 
minimum is zero when w, = 1 for all /, where w, is the weight of beamlet / of a 
plurality of radiation beams. 

33. The method of claim 30, wherein providing radiation comprises 
25 providing intensity-modulated radiation therapy. 

34. The method of claim 30, wherein providing radiation comprises 
providing Tomotherapy. 

30 35. A system for optimizing dose delivery of radiation comprising: 

computer processing means for determining an objective function to be 
used for mapping radiation beams to a body volume comprising at least one 



15 
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target volume, and at least one non-target volume, the objective function 
comprising a first term related to the at least one target volume and a second 
term related to the at least one non-target volume, the computer processing 
means being arranged to determine a minimum of the objective function 
5 whereby beams mapped so as to pass through the at least one non-target 
volume are limited such that the second term is zero only if the weights of 
beamlets passing through the at least one non-target volume are zero; and 

data communication means operably associated with the computer 
processing means for providing data to a radiation delivery apparatus for 
10 delivering radiation to the body volume based on the determined minimum of 
the objective function. 

36. The system of claim 35, wherein the dose delivery of radiation 
comprises intensity-modulated radiation therapy. 

15 

37. The system of claim 35, wherein the dose delivery of radiation 
comprises Tomotherapy. 

38. Computer readable storage having stored thereon computer program 
20 instructions executable on a computer system for causing the computer 

system to perform a method comprising: 

determining an objective function to be used for mapping radiation 
beams for a body volume comprising at least one target volume and at least 
one non-target volume, the objective function comprising a first term related to 

25 the at least one target volume and a second term related to the at least one 
non-target volume; and 

determining a minimum of the objective function whereby beams 
mapped so as to pass through at least one non-target volume are limited such 
that the second term is zero only if intensities if beamlets passing through the 

30 at least one non-target volume are zero. 
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39. The method of claim 10, wherein the importance parameter weighting 
each organ volume matrix is determined according to a function of position 
within the respective organ volume. 

40. The method of claim 10, wherein each importance parameter has a 
5 predetermined value. 
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